format long


% T = exp(T);

sigma2 = 0.66;
N2 = 3;
rho1 = 0.975;
[ZprobrouT, ZrouT] = rouwen(rho1, 0.0, sigma2, N2)
ZexprouT = exp(ZrouT)
ZprobrouT = transpose(ZprobrouT)
%sig_uncond=std(e)/sqrt(1-rho_Rouw^2);
sig_e = sigma2*sqrt(1-0.975^2)

Test = ZprobrouT;
F0 = zeros(1,N2);
F0(1,1) = 1;
tol = 1;
while tol > 10^(-8);
   F1 = F0*ZprobrouT;
   tol=max(abs(F1-F0));
   F0 = F1;
end;

%Normalize to 1:
ZrouTinv = F0./sum(F0);



rho = 0.8;
sig = 0.18;
N = 3;
mu = 0.0;
[ZT, ZprobT] = mytauchen(0.0, rho, sig, 3); %rouwen(rho,mu,sig,N);
ZexpT = exp(ZT)
ZprobT = ZprobT' % only for rouwen


F0 = zeros(1,N);
F0(1,1) = 1;
tol = 1;
while tol > 10^(-8);
   F1 = F0*ZprobT;
   tol=max(abs(F1-F0));
   F0 = F1;
end;

%Normalize to 1:
ZTinv = F0./sum(F0);


% % normalise to 1
% sumY = 0;
% for i = 1:N
%     sumY = sumY + ZTinv(i)*ZexpT(i);
% end
% for i = 1:N
%     ZexpT(i) = ZexpT(i)/sumY;
% end
% ZexpT






Mass = zeros(N*N2, 1);
wage = zeros(N*N2, 1);
Masstest = 0.0;
for i = 1:N
    for j = 1:N2
        Mass(i+(j-1)*N) = ZrouTinv(j)*ZTinv(i);
        Masstest =  Masstest + Mass(i+(j-1)*N);
        wage(i+(j-1)*N) = ZexprouT(j)*ZexpT(i);
    end
end
Masstest;

g = gini(Mass,wage)

std = var(ZexpT,ZTinv);

stop
% 
% rho = 0.87;
% sig = 0.295;
% N = 6;
% mu = 1.0;
% [ZT, ZprobT] = mytauchen(mu,rho,sig,N);


PP = [0.9737,0.0051,0.0212; 0.0524,0.9360,0.0116 ; 0.3,0.022,0.678];


F0 = zeros(1,3);
F0(1,1) = 1;
tol = 1;
while tol > 10^(-8);
   F1 = F0*PP;
   tol=max(abs(F1-F0));
   F0 = F1;
end;

%Normalize to 1:
PPinv = F0./sum(F0);



PPP = [0.9734, 0.0053, 0.0213; 0.0569, 0.9316, 0.0115; 0.4661, 0.0220, 0.5119];

PPP = [0.9756, 0.0052, 0.0192; 0.0518, 0.9398, 0.0084; 0.4305, 0.0239, 0.5456];

PPP = [0.9711, 0.0084, 0.0205; 0.0479, 0.9470, 0.0051; 0.4308, 0.0278, 0.5414];

F0 = zeros(1,3);
F0(1,1) = 1;
tol = 1;
while tol > 10^(-8);
   F1 = F0*PPP;
   tol=max(abs(F1-F0));
   F0 = F1;
end;

ZprobrouT



